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regions where superconductivity is accompanied by incommensurate magnetic order. 
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Introduction. The two-dimensional Hubbard model is a 
prototype system for competing order in layered transition 
metal oxide compounds. Shortly after the discovery of 
high temperature superconductivity in cuprates, it has 
been proposed as a key model describing the valence 
electrons in the copper-oxygen planes [T]. Indeed, the 
model exhibits the most prominent ordering phenomena 
observed in high-Tj, cuprates, namely antiferromagnetism 
and d-wave superconductivity [2]. 

While the magnetic order is captured already by 
conventional mean-field theory |3], superconductivity is 
fluctuation-driven and hence more subtle. Nevertheless, 
the emergence of d-wave superconductivity in the 2D Hub¬ 
bard model is nowadays well established [5]. In particular, 
unbiased evidence for superconductivity with sizable gaps 
at moderate interaction strengths has been obtained from 
functional renormalization group (fRG) calculations [IHS] , 
and from embedded quantum cluster methods at interme¬ 
diate and strong coupling [THIS]. 

The magnetic order is not necessarily of commensurate 
Neel type, that is, with antiparallel spin orientation be¬ 
tween adjacent sites. The possibility of magnetic order 
with generally incommensurate wave vectors away from 
(tt, tt) has been explored by several mean-field studies 
[HHH], and also by expansions in the limit of a small 
hole density, where fluctuation effects were taken into 
account |19I121) . Incommensurate magnetic order in the 
two-dimensional Hubbard model is also indicated by di¬ 
verging interactions and susceptibilities at incommensu¬ 
rate momenta in fRG flows However, the 

competition and possible coexistence of incommensurate 
magnetism and superconductivity has not yet been an¬ 
alyzed [23]. To do this, one needs to access the ordered 
phase in a framework that captures the fluctuations which 
generate d-wave superconductivity, allowing at the same 
time for a high momentum space resolution to distinguish 
the incommensurate ordering wave vector from (7r,7r). 
The latter requirement is an obstacle for cluster meth¬ 
ods, which have so far been restricted to commensurate 


antiferromagnetism jEHiniiia. 

In the fRG flow spontaneous symmetry breaking is 
signaled by diverging effective interactions at a critical 
energy scale Aj,. In principle, the flow can be continued 
below the instability scale Ac to compute order param¬ 
eters in the ordered phase [311251126]. However, this is 
rather complicated, especially in case of two or more order 
parameters. A drastic simplification occurs if the flow 
below the scale Ac is approximated by mean-field theory 
(MFT), since then order parameters can be computed 
without dealing with anomalous interactions [27]. While 
fluctuations above Ac are crucial for the generation of 
d-wave superconductivity, fluctuations below Ac are not 
expected to affect ground state order parameters signifi¬ 
cantly. Recently, a consistent fusion of the fRG flow above 
the scale Ac with MFT below Ac has been formulated 
and used to study the competition of antiferromagnetism 
and superconductivity in the Hubbard model |28j . In the 
purely superconducting regime of the ground state phase 
diagram, the pairing gap computed by fRG-|-MFT could 
be compared to previous results from a complete fRG flow 
with symmetry breaking (including fluctuations below Ac) 
[5]. The results obtained from the two methods agree very 
well. Away from half-filling, antiferromagnetic order was 
shown to be accompanied by microscopically coexisting 
d-wave superconductivity, in qualitative agreement with 
results from a previous fRG-|-MFT computation m and 
from cluster methods at stronger interactions [8Hini[i3]. 

Here we present novel results for the two-dimensional 
Hubbard model, allowing, for the first time, magnetic 
order with arbitrary wave vectors and superconductiv¬ 
ity. In addition to the previously established regions of 
Neel order coexisting with d-wave superconductivity, we 
find new coexistence regions where superconductivity is 
accompanied by incommensurate magnetic order. 

Method. Our results are based on the fRG-l-MFT 
method formulated in Ref. [28]. The calculation consists 
of three steps. First, the two-particle vertex is computed 
from a fRG flow integrated down to a scale Amf slightly 
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above the critical scale Ac. The flow is approximated 
by a one-loop truncation, with a static (frequency in¬ 
dependent) vertex parametrized via a decomposition in 
charge, magnetic, and pairing channels |23j with s-wave 
and d-wave form factors as described in Ref. [5]. The 
scale dependence defining the fRG flow is introduced by a 
soft frequency cutoff regularizing the bare propagator as 


Go(fco,k)= i sgn{ko) J-{ek-fi) 


( 1 ) 


We denote the spin and momentum dependence of the 
vertex at scale A by r^^_cr2;<T3,f74(ki,k2;k3,k4), where the 
indices 1,2 refer to outgoing, and 3,4 to incoming particles. 
Spin-singlet pairing is driven by the singlet component of 
the vertex in the Cooper channel (zero total momentum), 

Vick' = ^rA“-(k,-k;-k',k'), ( 2 ) 

where T^ = Magnetic instabili¬ 

ties with a momentum Q are generated by the magnetic 
component of the vertex with momentum transfer Q, that 
is, 

c/Q;kk' = r^,r^;.,_.(k -I- Q, k'; k' + Q, k) . (3) 

Note that T^ due to spin 

rotation invariance. 

In the second step, the irreducible vertices in the rele¬ 
vant channels are computed by solving the corresponding 
Bethe-Salpeter equations [28l |29] . The irreducible pairing 
vertex I4;k' is obtained from the equation 


Vick' = Vick' 

- f I4pG^“"(po,p)G^“M-po,-p)Vpk', (4) 
Jp 

where G^“'' {po, p) is the propagator at the scale Amf, and 
is a short-hand notation for J ^. The irreducible 
magnetic vertex C/Q;kk' is determined from 


U, 


Q;kk' 


= G, 


Q;kk' 


+ f GQ,kpG^“-(po,p)G^“"(Po,P+Q)GQ: 

^ P 


pk' 


(5) 


In this work we neglect self-energy contributions to the 
flow, so that G^“'' is just the regularized bare propagator 

G Amf 
0 

Finally, the superconducting and magnetic order param¬ 
eters are computed by solving the mean-field Hamiltonian 
with the irreducible vertex parts Vick' and GQ;kk' as effec¬ 
tive interactions. The superconducting order parameter 
is given by the gap function 

*l< = /(§)Uk.(pcO, (6) 


where pk = O-k^Okt is the Cooper pair annihilation op¬ 
erator. The phase of the superconducting order can be 
chosen such that Ak is real. Concerning the magnetic 
order, we restrict our analysis to spiral states, as described 
by the order parameter 

^k = GQ;kk'(?Tlk') , (7) 

where mk = a|^.^ak-i-Q^. For Q = (tt, tt) the spiral order is 
equivalent to Neel order with a staggered magnetisation in 
the a:?/-plane. There are several mean-held studies of spiral 
magnetic order in the two-dimensional Hubbard model 
miiiHi- A frequently discussed alternative is collinear 
order [I1II5], especially in combination with pronounced 
charge stripes In a spiral state the amplitude of 

the magnetization is homogeneous, only the orientation 
varies. Hence, the magnetic order entails an energy gain 
everywhere in the system. By contrast, collinear order 
necessarily involves regions with a reduced magnetization, 
where the energy gain from the order is also reduced. 
Incommensurate collinear order is therefore expected to 
be favorable only in the form of sharply dehned charged 
domain walls between antiferromagnetic regions [HIH]. 
In the moderate interaction regime investigated here such 
pronounced real space prohles seem unfavorable since they 
cost a lot of kinetic (hopping) energy. 

A mean-held decoupling of the reduced effective inter¬ 
actions yields the mean-held Hamiltonian 

FImF — 'y ^ f'kaj^^Ukcr 

k,0' 

+ ^ ^k^k + ^k^k - (^k) 

k 

+ ^Akp|, +A^Pk-A^(;?k)- (8) 

k 

For the Hubbard model with nearest and next-to-nearest 
neighbor hopping on a square lattice, the dispersion 
relation is Ck = —2t{coskx + cosky) — At' cos kx cos ky. 

Using Nambu spinors d^k = (^k^)®-ki>“k+Q;>“-k-Qt)> 
the mean-held Hamiltonian can be written in the form 
Hmf = 5 X^k 'I']{jA4k4'k+const, where Aik is a hermitian 
4x4 matrix. Hmf can thus be diagonalized by a 4 x 4 
unitary (generalized Bogoliubov) transformation, and the 
resulting gap equations can be solved numerically by iter¬ 
ation. Occasionally two distinct locally stable solutions 
of the gap equations are found. One then has to compute 
the corresponding free energies to discriminate globally 
stable from metastable states. 

Results. We have computed the magnetic and super¬ 
conducting order parameters in the ground state of the 
two-dimensional Hubbard model at and near half-hlling 
for weak to moderate interaction strenghts G/I = 2, 3,4 
for a small next-to-nearest neighbor hopping amplitude 
t'/t = —0.15, and for U/t = 3 also for the special particle- 
hole symmetric case t' = 0. We now present our results 
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Figure 1: (Color online) Amplitudes of magnetic and super¬ 
conducting gap functions in the ground state of the two- 
dimensional Hubbard model as a function of the electron 
density n, for t/ = 3 and t' = —0.15. The amplitude of the 
magnetic gap function obtained from a purely magnetic so¬ 
lution is also shown for comparison, with corresponding hole 
and electron pockets at n = 0.96 and n = 1.05 in the left and 
right insets, respectively. First order transitions (see broken 
lines near n = 0.9 and n = 1.06) lead to small density intervals 
where no homogeneous solution exists. 


for f//t = 3 and t' jt = —0.15, and subsequently discuss 
similarities and differences found for the other choices of 
U and t' m- In the following we set t = 1. 

In Fig. 1 we show the amplitudes of the magnetic and 
superconducting order parameters, A = maxk Ak and 
A = maxk Ak, as a function of the electron density. The 
magnetic order is of Neel-type (Q = (tt, tt)) at half-filling. 
The Neel order persists on the electron-doped side (n > 1), 
and also for small and moderate hole doping (0.9 < n < 1). 
The stability of commensurate antiferromagnetic order 
on the electron-doped side is generally expected [55]. For 
n < 1 excitations in hole-pockets near the Brillouin zone 
diagonal could destabilize the commensurate state even 
for small hole doping for small |t'| and/or large U [20] . 
At n — 0.9 a first order transition to an incommensu¬ 
rate spiral state occurs, with a wave vector of the form 
Q = (tt — 27r77,7r), or equivalent (by symmetry) wave vec¬ 
tors (tt + 27r77, tt) or (tt, tt ± 2TTr]). The magnetic order is 
completely suppressed by superconductivity at van Hove 
filling (n = 0.87), but then reemerges for lower densities. 
The magnetic transition at n = 0.73 is of weak first-order 
type. The incommensurability ry is plotted as a function of 
the density in Fig. 2. It jumps from zero to a small finite 
value at the commensurate-incommensurate transition, 
and then increases monotonically upon further doping 
until the magnetic order disappears at n = 0.73. The 
chemical potential /J.(n), which is also plotted in Fig. 2, ex¬ 
hibits a discontinuity due to the charge gap at half-filling. 
In other words, the density n is pinned to half-filling 
for a range of chemical potentials between —0.37 and 
—0.245. Comparing to the purely magnetic solution (ex¬ 
cluding superconductivity), which is also shown in Fig. 1, 
one can see that superconductivity has little influence 



Figure 2: (Color online) Incommensurability rj and chemical 
potential p as function of the density for the same parameters 
as in Fig. 1. 


on the magnetic order in the regime around half-filling 
where it is commensurate. The incommensurate magnetic 
order on the hole-doped side is strongly suppressed by 
pairing in the vicinity of van Hove filling. At van Hove 
filling, superconductivity eliminates the magnetic order 
completely. By contrast, in the overdoped regime well 
below van Hove hlling, incommensurate magnetic order 
and superconductivity coexist without suppressing each 
other significantly. 

The pairing gap Ak is finite for all densities except at 
half-filling, where the Fermi surface is fully gapped by the 
antiferromagnetic order. Hence, away from half-filling, 
magnetic order always allows for coexisting supercon¬ 
ductivity. This is easily understood as follows. Doping 
the half-filled antiferromagnet by additional electrons or 
holes leads to electron or hole pockets (see the insets of 
Fig. 1). The ubiquitous attractive d-wave pairing inter¬ 
action inevitably generates a Cooper instability at the 
(small) Fermi surfaces enclosing these pockets, and thus 
superconductivity. The hole pockets for n < 1 are cen¬ 
tered around the nodes of the pairing gap Ak, while the 
electron pockets for n > 1 are in the antinodal region 
near (tt, 0) and (0,7r), where the gap is maximal. Hence, 
the onset of pairing around half-filling is much steeper on 
the electron doped side (see Fig. 1). This is in agreement 
with recent spin-fluctuation exchange calculations in the 
weak coupling regime |39j . but differs from the behavior 
found in a strong coupling analysis of electrons moving 
in an antiferromagnetic spin background. In the latter 
case, pairing is mediated mostly by transverse spin fluctu¬ 
ations (magnons), which couple very weakly to electrons 
in the antinodal region, so that the pairing interaction 
is very small for low electron doping m- By contrast, 
at moderate coupling also longitudinal spin fluctuations 
contribute, and yield a sizable pairing interaction in the 
antinodal region [55] . 

The maximal gap amplitude at ’’optimal” hole-doping is 
significantly larger than the maximal gap on the electron 
doped side, in agreement with the gap hierarchy in cuprate 
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-- E(A,A)-E(0,0) -B- E(A,A)-E(0,A') 

Figure 3: (Color online) Total condensation energy E(A, A) — 
i?(0,0) and relative magnetic condensation energy E(A,A) — 
E{0, A') as a function of the chemical potential /r. The inset 
is a zoom into the region with incommsurate magnetic order. 
The model parameters are the same as in Fig. 1. 


superconductors. The maximal gap on the hole-doped 
side is situated slightly above van Hove filling, where the 
magnetic order is already quite weak. The gap decreases 
smoothly in the ’’overdoped” regime, due to a decrease 
of the pairing interaction and the density of states. 

The leading order parameter is often guessed from the 
leading divergence of the effective interactions in the fRG 
flow upon approaching the critical scale Ac [4]. This is 
usually correct, but there are exceptions. In particular, 
around van Hove filling, we find superconductivity as the 
dominant order, although the leading divergence occurs 
actually in the magnetic channel (see Ref. [S] where the 
fRG flow was studied for the same parameters). 

A major novel result of our work is the coexistence of 
superconductivity with incommensurate magnetic order 
in a regime where a purely superconducting state was 
expected, since the hole-doping is already quite large. 
It is therefore interesting to look at the condensation 
energy gained by symmetry breaking in the incommen¬ 
surate regime. In Fig. 3 we plot the total condensation 
energy E(A, A) — E{0, 0) gained by magnetism and pair¬ 
ing. We also plot the additional energy gain due to the 
magnetic order compared to a purely superconducting 
state, E{A, A) — E{0, A'), where A' is the pairing gap in 
the absence of magnetism. One can see that this addi¬ 
tional energy gain is tiny in the incommensurate regime 

< —0.57), even for densities where the size of the mag¬ 
netic order parameter is comparable to the pairing gap. 
Hence, the magnetic order in this regime is extremely 
delicate, a kind of ’’gossamer magnetism”, reminiscent of 
Laughlin’s ’’gossamer superconductivity” in lightly doped 
Mott insulators m- The quasi-degeneracy of a purely 
superconducting state and a state with coexisting mag¬ 
netism can be expected to lead to intriguing fluctuation 
effects. 

Let us now compare the above results to those obtained 


for other choices of I/, as presented in our Supplementary 
Material [37] ■ For U = 2 and t' = —0.15, the ground state 
is purely superconducting for all densities, with a small 
d-wave gap which is maximal slightly above van Hove 
filling. Magnetic order occurs only if superconductivity is 
switched off. This is expected in the weak coupling limit 
in the absence of nesting (for finite t'). For [7 = 4 and 
t' = —0.15 there is no homogeneous solution in the density 
range 0.88 < n < 1, which includes also van Hove filling. 
Hence, a system with an average density in that interval 
will undergo phase separation in regions with distinct den¬ 
sities n = 1 and n = 0.88, or form a more complex type 
of order. For smaller densities incommensurate magnetic 
order coexisting with d-wave superconductivity is found, 
again with a tiny energy gain from the magnetic order. 
Hence, our main result, the coexistence of superconductiv¬ 
ity with a very delicate incommensurate magnetic order 
at sizable hole-doping, is robust with respect to an in¬ 
crease of U. On the electron-doped side, there are no 
qualitative differences compared to 1/ = 3, except for a 
tiny incommensurate region at the edge of the magnetic 
regime. 

For t' = 0, the Hubbard model is particle-hole symmet¬ 
ric, that is, the properties for densities n and 2 — n are 
equivalent. Due to perfect nesting the system is a Neel 
antiferromagnet at half-filling for any U > 0. There is 
no homogeneous solution in a density range around half¬ 
filling, both on the electron and hole-doped side. For larger 
doping there is a small region exhibiting incommensurate 
magnetism in coexistence with d-wave superconductivity, 
and a broader purely superconducting region. While the 
magnetic order parameter at half-filling has almost the 
same size for t' = 0 and t' = —0.15, the largest achievable 
pairing gap (at optimal doping) is much smaller for t' = 0. 
A sizable next-to-nearest neighbor hopping thus helps to 
promote superconductivity with a large gap. This was 
already revealed in a previous fRG study |S] and in a re¬ 
cent quantum cluster calculation na, and is in qualitative 
agreement with the empirical trend in cuprates |42j . 

Summary. We have analyzed the competition between 
magnetism and superconductivity in the ground state of 
the two-dimensional Hubbard model, including the pos¬ 
sibility of incommensurate spiral magnetic order. Using 
a combination of fRG and mean-field theory, fluctuation 
driven order is captured without any bias for a specific 
instability. Gharge, spin, and pairing channels are treated 
on equal footing. Away from half-filling, magnetic order 
always coexists with superconductivity, as a consequence 
of a Cooper instability in electron or hole pockets. For 
E < 0, both magnetism and superconductivity exhibit a 
pronounced particle-hole asymmetry. On the hole-doped 
side, superconductivity has larger maximal gaps and it 
coexists with incommensurate magnetism at moderate 
doping, except at van Hove filling, where magnetic or¬ 
der is fully suppressed by pairing. The incommensurate 
magnetic order is “gossamer-like” in the sense that it is 
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stabilized only by a tiny energy gain with respect to a 
purely superconducting state. Rather fragile incommen¬ 
surate magnetic order has actually been observed at the 
bottom of the superconducting dome in La 2 _xSrxCu 04 
[33] and YBa2Cu306-i-x [HI US]- Suppressing supercon¬ 
ductivity by a strong magnetic field would stabilize that 
order. Indeed, in a very recent high field experiment a Hall 
coefficient consistent with hole-pockets arising from pos¬ 
sible spin-density waves was observed in YBa 2 Cu 306 +x 
at much higher doping than previously [46] . 

Our analysis was restricted to weak and moderate inter¬ 
action strengths. The method may be extended to strong 
interactions by using dynamical mean-field theory as a 
starting point for the fRG flow m- 
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COEXISTENCE OF INCOMMENSURATE MAGNETISM AND SUPERCONDUCTIVITY IN THE 
TWO-DIMENSIONAL HUBBARD MODEL: SUPPLEMENTAL MATERIAL 

H. Yamase, A. Eberlein, and W. Metzner 


In this Supplemental Material we provide details on the computation of the irreducible vertices, and we present 
additional results for alternative choices of the model parameters. 


Computation of irreducible vertices 

Here we discuss some technical aspects concerning the computation of the irreducible vertices Vkk' and C/Q;kk' from 
Vkk' and C/Q:kk', respectively. The k and k' dependences are discretized by partitioning the momentum space in 
patches, that is. 


Hck' — ) 

C^Q;kk' = Ejj'C^Q;H'0i(k)0i'(k'), 


where 0j(k) is the characteristic function of the patch labelled by j, that is, 0j (k) = 1 if k is on patch j and 
0j(k) = 0 otherwise. The solution of the B et he-S alp ether equations is thus reduced to a matrix inversion. The 
matrices Vjj> and C/q^j' are real and symmetric in j and j'. Hence, they can be decomposed into eigenmodes as 
Vjj' = ^ with real eigenvalues and real normalized eigenvectors. The 

leading instabilities correspond to the most negative eigenvalues. However, the discretized irreducible vertices Vjji 
and Uq-jji obtained from the Bethe-Salpeter equations contain eigenmodes with large negative eigenvalues which 
correspond to eigenmodes with large positive eigenvalues of the full vertices. To remove these unphysical instabilities, 
we project the vertices on the subspace with negative eigenvalues. 


TJI - V' 


(S2) 


where the primed sums are restricted to negative eigenvalues. The irreducible vertices V'j, and Uq.jj, computed from 
the projected full vertices do not contain any unphysical instabilities, and the leading physical instabilities are not 
affected by the projection. 


Results for alternative choices of model parameters 

In this section we present results for alternative choices of the model parameters, which complement the results for 
t'/t = —0.15 and U/t = i presented in the main paper. In the following we set t = 1. 

In Fig. SI we show the amplitudes of magnetic and superconducting order parameters, A = maxkAk and A = 
maxk Ak, as a function of the electron density for t' = —0.15 and two distinct interaction strengths, {7 = 2 (left) and 
U = 4: (right). For {7 = 2, the minimal energy is reached by a purely superconducting solution at all densities. The 
d-wave gap function has a sizable amplitude only in a density range between n = 0.83 and n = 0.96, with a maximum 
slightly above van Hove filling (rivH = 0.875). In the absence of superconductivity (setting Ak = 0 in the mean-field 
equations), a magnetic solution emerges in a small density range around van Hove filling. The magnetic order is 
commensurate for densities n > 0.879 and incommensurate below. The commensurate-incommensurate transition 
is thus situated slightly above van Hove filling. It is discontinuous, with a very small jump of the magnetic order 
parameter (not visible in the figure). The deviation of the ordering wave vector from (tTjTt) is quite small (ry < 0.04). 
The incommensurate magnetic transition near n = 0.86 is continuous, the commensurate transition at n = 0.895 is 
discontinuous. 

For {7 = 4 the model exhibits Neel order with a sizable magnetic order parameter at half-filling. The Fermi surface is 
fully gapped, and the corresponding charge gap leads to a pinning of the density at half-filling for a range of chemical 
potentials from /i = —0.64 to /r = —0.036. While the ground state should be invariant in that range, our approximate 
results for the magnetic order parameter Ak depend to some extent on the choice of p,. In Fig. SI we have plotted 
A as obtained at p = —0.036. At the other edge, A is about 10 percent larger. This rather strong increase occurs 
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Figure SI: Amplitudes of magnetic and superconducting gap functions in the ground state of the two-dimensional Hubbard 
model for a next-to-nearest neighbor hopping amplitude t' = —0.15 and interaction strengths U = 2 (left) and ?7 = 4 (right). 
The amplitude of the magnetic gap function obtained from a purely magnetic solution is also shown for comparison. 




near the van Hove singularity of the non-interacting system at /t = —0.6. For C/ = 3 the artificial dependence on /i 
is much smaller (below one percent). At /r = —0.64 there is a first order transition between the Neel state and a 
phase with incommensurate magnetic order coexisting with d-wave superconductivity. The transition is accompanied 
by a large density drop to n = 0.88. No homogeneous solution exists for densities between n = 0.88 and n = 1, so 
that systems with densities in that range must undergo phase separation. The magnetic order parameter vanishes 
continuously at n = 0.66. In the entire density range 0.66 < n < 0.88 the incommensurate magnetic order coexists 
with d-waye superconductivity. While the amplitudes of the two order parameters are generally comparable in this 
coexistence region, the energy gain resulting from the magnetic order is tiny compared to the gain obtained by pairing 
(see the inset of Fig. S2). This “gossamer-type” magnetism was already observed for U = 3 (see Fig. 3). On the 



Figure S2: Total condensation energy E{A, A) — E{0, 0) and relative magnetic condensation energy E{A, A) — E{0, A') as a 
function of the chemical potential fi, for t' = —0.15 and [7 = 4. 


electron-doped side the commensurate magnetic order persists in a wide density range. There is a first order transition 
to an incommensurate state with a tiny incommensurability at n = 1.15, near the edge of the magnetic region. The 
magnetic order parameter vanishes continuously at a slightly higher density. The magnetic order on the electron-doped 
side also coexists with d-wave superconductivity, with a steep onset of pairing near half-filling. The pairing gap remains 
sizable in a small density range outside the magnetic regime. 

We now turn to the special case with pure nearest-neighbor hopping (f = 0) for comparison. The Fermi surface at 
half-filling is perfectly nested in that case. Hence, at half-filling the system is a Neel antiferromagnet for any 17 > 0. 
The system is particle-hole symmetric, that is, the properties for densities n and 2 — n are equivalent, so that we can 
restrict ourselves to the hole-doped region. We have computed the magnetic and superconducting order parameters 
for U = 3. At half-filling the system is purely antiferromagnetic with a magnetic order parameter A = 0.36. No 
homogeneous solution is found in the density range 0.875 < n < 1. Hence, a system with an average density in that 
interval will undergo phase separation in regions with distinct densities n = 1 and n = 0.875. The magnetic and 
superconducting order parameters for n < 0.875 are shown in Fig. S3. For n < 0.875 there is a small region exhibiting 
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incommensurate magnetism in coexistence with d-wave superconductivity, and a broader purely superconducting 
region below. The incommensurability is small (77 = 0.041) and almost density-independent in the narrow regime 
where the incommensurate magnetic order exists. 



Figure S3: Amplitudes of magnetic and superconducting gap functions in the ground state of the two-dimensional Hubbard 
model with t' = 0 and t/ = 3 in the density range 0.848 < n < 0.875. The amplitude of the magnetic gap function obtained 
from a purely magnetic solution is also shown for comparison. 


The extent of phase separated regions obviously depends sensitively on the model parameters. The issue of phase 
separation in the Hubbard model is discussed controversially in the literature. Mean-field theory yields phase separation 
in the ground state at = 0 near half-filling for any U > 0, while a finite t' is found to stabilize homogeneous solutions 
at weak interactions both in the electron- and hole-doped regime [SI]. This is in full agreement with our results. At 
strong coupling, phase separation at = 0 is found by a dual-fermion approach [S2], and on the hole-doped side 
at t' = —0.3 within a variational cluster approximation [S3]. By contrast, in a dynamical cluster calculation phase 
separation was found only for a positive t' on the hole-doped side (equivalent to a negative t' and electron-doping) [S4]. 
Quantum Monte Carlo calculations and exact diagonalization do not provide evidence for phase separation at all, but 
this may well be due to finite size effects. 
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